Formal Solution of a Class of Reaction-Diffusion Models: Reduction to a 

Single-Particle Problem 
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We consider the trapping reaction A + B ^ B m space dimension d < 2. By formally eliminating 
the B particles from the problem we derive an effective dynamics for the A particles from which the 
survival probability of a given A particle and the statistics of its spatial fluctuations can be calculated 
in a rather general way. The method can be extended to the study of annihilation/coalescence 
reactions, _B + -B ^ or B, in ci = 2. 
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First-passage problems involving more than a few de- 
grees of freedom are notoriously difficult to solve |j. 
In this Letter we introduce a technique that enables one 
to solve a class of first-passage problems involving an in- 
finite number of degrees of freedom. For definiteness we 
develop the method in the context of the "trapping reac- 
tion", A + B B, but the applications are more general, 
as emphasized in the latter part of the Letter. The main 
result of our approach is to reduce the problem to one 
described by a single degree of freedom whose late-time 
behavior can be extracted analytically. 

The asymptotic dynamics of the trapping reaction has 
been a long-standing puzzle. The main question is how 
the density of A particles decreases with time. A re- 
lated problem, much studied in the context of chemical 
kinetics [||, Q is the two-species annihilation reaction, 
A + B ^ 0, with initial densities pa{0) < Pb{0). This is 
equivalent to the trapping reaction at late times, when 
PA{t) ^ Psit) and pB{t) is essentially constant. Again, 
the standard problem is to compute the asymptotic form 
of the A-particle density pAit) or, equivalently, the prob- 
ability Q(t) that a given A particle survives until time t. 
Since the particles do not interact with other particles of 
the same species, to compute Q{t) it suffices to consider 
a single A particle moving in an infinite sea of B particles 
with density p {— pb)- 

Traditional approaches to this type of problem treat 
the A particle as an absorbing boundary for the B par- 
ticles. Unfortunately, for an arbitrary A-particle trajec- 
tory this absorbing-boundary problem cannot be solved. 
In this Letter we introduce a different approach in which 
we treat the A and B particles as if they were non- 
interacting. We exploit the initial condition that each 
B particle is randomly located anywhere in the system 
to show that certain 'events', where a B particle meets 
the A particle for the first time (remember that we are 
treating them as non-interacting, so they can meet more 
than once) have a Poisson distribution, i.e. the proba- 
bility pn that n such events have occurred up to time t 
is given by p„ = (/^"/n!) exp(— /^), where the mean, p,, 
of the distribution is a functional^ p[z\ of the trajectory 
z{t), < r < of the A particle. The probability that 



the trajectory z(t) has survived, in the original interact- 
ing problem, is simply po[z\ = cxp(— /i[z]). Finally, Q{t) 
is obtained by averaging exp(— /i[z|) over all possible A- 
particle trajectories z(t) with the appropriate (Wiener) 
measure, exp [-{II AD') J^dT{dz/dT)% where D' is the 
A particle's diffusion constant. In this way, the B parti- 
cles have been eliminated from the problem, and one has 
an effective A-particle dynamics described by the Wiener 
measure and the functional p[z\. The final step, that 
makes further analytical progress possible, is the obser- 
vation that the path integral over z(t) is dominated at 
late times by a single A-particle trajectory. 

The main results of this approach are: (i) the tra- 
jectory where the A particle is stationary is proved to 
be the dominant trajectory and determines the asymp- 
totic form of the A particle's survival probability 
Q{t) ~ exp(— Adt''^^) for d < 2 (with a logarithmic cor- 
rection in d ~ 2), where is a calculable constant 
and d is the dimensionality of space; (ii) typical fluc- 
tuations of surviving A-particle trajectories around this 
dominant path have variance {z^{t)) ^ t(2-d)/2 ^^j. 2 
(iii) exact results are obtained for Q{t) and the form of 
the dominant path in a system with a non-uniform ini- 
tial density of B particles; (iv) this approach provides 
a powerful method for calculating first-passage proper- 
ties for a deterministically moving boundary z{t) ; (v) the 
method also provides a formalism for calculating Q{t) 
in the highly nontrivial situation where the B particles 
themselves interact, e.g. i? -|- B ^ 0, at least in d = 2 
where the density correlations induced by these reactions 
are negligible. 

We begin by deriving the Poisson property that plays 
a central role in the analysis. We consider a finite vol- 
ume, V , containing N — pV _B-particles (diffusion con- 
stant D), randomly distributed within it, and a single 
A particle (diffusion constant D'), initially located at 
the origin. Let z{t) be the A particle's trajectory, and 
let P{x, t) be the probability that a given B particle, 
starting at x, has met the A particle before time t. The 
average of this quantity over the initial position, x, is 
{1/V) JydV P{x,t) = R{t)/V, where R{t) is an implicit 
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functional of z{t). The probability that n distinct B 
particles have met the A particle, averaged over their 
initial positions, is p„(t) = (^) - R/V)^-". 
Taking the limit ^ oo, ^ ^ oo, with p = N/V 
and n held fixed, one recovers the Poisson distribution, 
Pn — (/^"/'^O 6xp(— /^), with fi = pR. 

One can derive an equation for the functional p[z\ by 
calculating, in two ways, the probability density to find a 
B particle at the point z{t) at time t. First, since the par- 
ticles are treated as non-interacting, and the B particles 
start in a steady-state configuration of uniform density p, 
this probability density is just p. Secondly, from the Pois- 
son property, the probability that a B particle (i.e. any B 
particle) meets the A particle for the first time in the time 
interval [t' ,t' + dt') is p{t')dt' . The probability density 
for such a particle to subsequently arrive at z{t) at time 
t is given by the diffusion propagator G{z{t), t\z{t'),t') = 
\\-KD{t-t')\-'^l'^ exp{-[z(t)-z(i')]V4-D(i-^')}- Equat- 
ing the results from these two methods gives our funda- 
mental equation. 
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p^ \ dt' fi{t')G{z{t),t\z{t'),t') 
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which is an implicit equation for the functional p[z\ (not- 
ing that p{t — 0) = 0, since no B particle can meet the 
A particle in zero time). Finally, Q{t) = {exp{— p[z\)) ^ , 
where the average is over all paths z{t) weighted with the 
Wiener measure. 

As a first application of this equation we prove that 
the trajectory z = is the dominant path, i.e. that it 
gives the smallest possible value of p[z\ for all t. This 
fimction, po{t), satisfies Eq. (0) with z = 0: 



dt' pQ{t')[ATiD{t-t')]-'^/^ 
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wheTeK{ti,t2) = l-exp{-[z(ti)-£(i2)]V4-D(ii-i2)}- 
Eq. (|^) is 'implicit' because the full p appears on the 
right-hand side. Now note that K{ti,t2) > and p>0 
(because p{t) is the mean number of different B parti- 
cles that have met the A particle up to time t - clearly 
a non-decreasing function). Therefore pi[z\ > for all 
paths z{t), with equality when z{t) — for all t. It fol- 
lows that Q{t) = (exp(— /io — pi))z < exp[— /io(t)]. This 
rigorous upper bound for Q{t), combined with the iden- 
tical rigorous lower bound derived in 0, proves that the 
asymptotic form of Q{t) is the same as for the 'target' 
problem, where the A particle is stationary, for all d < 2. 
The interpretation of this result is that, since p is large 
for f ^ oo {p ^ f^/^), the path-integral for Q{t) is domi- 
nated by the path that minimizes p, i.e. we are essentially 
evaluating the path integral by the method of steepest 
descents. Small fluctuations around the dominant path 
will determine the corrections to the asymptotic form. 

We next compute the probability distribution, P{z,t), 
of the position z of the A particle at time t, given that 
it survives. Numerical studies |^ suggest that, in d = 
1, (^^(t))^/^ ~ f^, with (/) = 0.25 - 0.3, while similar 
studies in d = 2 are inconclusive. Our methods give 
</> = 1/4 in d = 1 and = (2 - d)/4 for all d < 2. The 
technique is to expand pi [z\ , given by Eq. (^) , to order 
to compute the variance of the Gaussian fluctuations 
around the dominant trajectory z{t) — 0. To this order 
one can replace p{t2) on the right-hand side by po(t2) — 
X^t^"^ ^ and expand the function K(ti,t2) to order z'^. 
Specializing to d = 1, the result is, at time t, 



By inspection, po{t) must have the form po{t) — Xdf^^^ 
(for d < 2) , in order that the right-hand side be indepen- 
dent of t. Substituting this form in (^), and evaluating 
the integral, gives 
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d< 2, 



(3) 



while for d = 2 one finds, for t oo, pa{t) ATrpDt/ hit 
1^ . The corresponding A-particle survival probability is 
Qo{t) — exp[— ^o(^)]- This simple case of a static A parti- 
cle is sometimes called the 'target annihilation problem', 
and our method reproduces the known results for that 
problem |^ in a very simple way. To prove that z(t) = 
gives the global minimum of p[z\ we write p = po + Pi 
in (|l|). This equation can then be rearranged, with the 
help of Laplace transform techniques, to give an implicit 
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The probability distribution for z at time t is given, up 
to an overall normalization, by the path integral 



P(z,<) = I Pz(i)cxp (-^^ 



drP^r) — pi [z 



(6) 

where the integral is over all paths satisfying the bound- 
ary conditions z(0) = 0, z{t) — z. 

The path integral has the form / I?z(t) exp(— 5'[z]), 
where the 'action' S[z\ = St[x\ + Sv[z\ is a quadratic 
functional of z{t), where St[z\ = (1 / 4:D') dr P (t) , 
and Sv[z] = Since the integrand is Gaussian, 

the z-dependence of the path integral is exactly cap- 
tured by the extremal path connecting z(0) = and 
z{t) = z. Although we have been unable to find this 
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path analytically, power counting on the two terms sug- 
gests Sv ~ Xxz^ ID\ft and St ~ z'^/D't, so that Sy 
dominates at large t. A more careful analysis ^ confirms 
the dominance of Sy and the form Sy = Cy^iz^ / Dy/i, 
where cy is a pure number. The probability weight for 
the fluctuations z{t) of surviving trajectories is therefore 
Gaussian, with variance {z'^it)) = {D/2cvXi)Viior large 
t, with the value of cv determined by the extremal path. 
Note that this result is independent of the A-particle dif- 
fusion constant D' . If D' = 0, of course, there are no 
fluctuations, but for any D' > the variance becomes 
independent of D' at large enough t. Also, since typically 
z ^ t^/^, the next term in the expansion of the function 
K{ti,t2) in is of relative order z'^ / Dt ~ i"^/^ so it 
is negligible at large t. The distribution of z is therefore 
exactly Gaussian, at least in the 'scaling limit' z — > oo, 
t oo, with z/t^/^ fixed. Similar arguments |^ give the 
generalization {z'^(t)) - (D/A<j)<(^"'*'/^ for d<2. 

We turn now to a related problem with a non-trivial 
dominant path that can be exactly determined. Con- 
sider, in d = 1, a system where the density of B parti- 
cles at t = has different values, pL and pR, to the left 
and right of the A particle. The derivation of an equa- 
tion for p,[z\ proceeds exactly as before, except that the 
probability density to find a B particle at the point z at 
time t in the noninteracting system, which appears on 
the left of Eq. (|l|), has to be recalculated. In terms of the 
diffusion propagator G introduced earlier, this probabil- 
ity is (quite generally) Pb(z, t) = dx p{x) G{z, t\x, 0) 
where p{x) is the initial _B-particle density at position x. 
When p{x) — p, a. constant, one finds Psiz^t) — p, as 
before. When p{x) — pL for x <{) and pR for a; > 0, the 
generalized version of (Q) becomes 



P 



1 - A erf 



z{t) 



dt' p{t')G{z{t),t\z{t'),t') 



(7) 

where p = (pl + Pr)/'^ is now the mean density, A = 
{pl — Pr)I{pl + Pr) is a measure of the left-right asym- 
metry, and erf(x) is the error function. 

Physical intuition suggests that, because of the asym- 
metry, surviving A-particle trajectories will tend to be 
those that drift into the region (the right, say) where the 
-B-particle density is initially smaller. Upper and lower 
bounds have been derived earlier |ic| ] for the asymptotics 
of the A-particle survival probability, Q{t), which show 
that it has the asymptotic form Q{t) ~ exp[— (7(A)Ai\/t], 
where f;(0) = 1 for consistency with the symmetric case, 
PL — PR- This form for Q{t) shows that p[z\ for the opti- 
mal path has the time-dependence p oc VI. Both sides of 
(^ can then be rendered time-independent by the choice 
z{t) = aV ADt for all t < t, where a is a constant to be 
determined. A more detailed analysis Q shows that the 
dominant path is indeed of this form. 

Putting p{t) = g{A)XiVi and z(t) = aV^Dr, in (0), 




ln(p Dt) 5 



10 



FIG. 1: Time-dependence of g{t) = —\nQ{t)/X\y/t for den- 
sity asymmetries A — 1/3, 3/5 and 7/9 (top to bottom). The 
horizontal lines show the asymptotic value of g in each case, 
while the inset shows this value as a function of A. 



and evaluating the integral on the right-hand side gives 



S(A) 



exp(— a^) 
1 - erf 2 (a) 



[1 - A erf (a)]. 



(8) 



The final step is to minimize the right-hand side with 
respect to a to obtain the optimal path. This can be 
done numerically. The resulting g{A) is shown for A > 
as the inset in Figure 1 (note that, by symmetry, g{A) is 
symmetric around A = 0). It clearly satisfies the bounds 
1 - |A| < g{A) < 1 derived in Q. 

In the main part of Fig. 1, numerical results for g{t) — 
— \nQ(t)/Xiy/i, obtained using the algorithm of ref. 
are displayed, with D ~ 1/2 and pL — 0.5 in all cases, 
while PR = 0.25,0.125,0.0625 (top to bottom on the 
right), corresponding to A = 1/3,3/5,7/9 respectively. 
The horizontal lines on the right show the asymptotic 
values obtained from the inset for the corresponding val- 
ues of A. The slow approach to asymptopia is similar to 
that observed ||, 10 1 in the symmetric case (A = 0). 

As a bonus, the same calculation solves the first- 
passage problem of B particles with density pl for cc < 
and PR for a; > moving in the presence of a de- 
terministically moving absorbing boundary located at 
x(t) = aV^Dt. The probability that no particle has 
reached the boundary up to time t has the simple form 
Q{t) = exp[-c/(A)Aiv^], with 5(A) given by (|). We 
are not aware of any other way of obtaining this result. 
Extensions to deterministically moving absorbing bound- 
aries in dimension d > 1 are also possible Q . 

In the final part of this Letter, we apply this ap- 
proach to a nontrivial problem with d ~ 2. Consider 
the annihilation/coalescence reaction B + B with 
probability l/{q — 1) and B + B ^ B with probability 
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{q — 2)/{q — 1). The density of B particles is known to 
decay as p{t) = ad[{q— I) / q]{Dt)~'^/^ for d <2, where 
is a universal constant ||ll|, |l^ equal to l/27re for d ^ 2, 
with e = 2 — d. For d — 2, a logarithmic correction is 
obtained, p{t) ~ \nt/t. Now suppose that one of the 
B particles is tagged and relabeled an A particle, with 
diffusion constant D' . We consider the probability (the 
'walker persistence' probability Q{t), that the A 

particle has not met any B particle up to time t. The 
limit d —^ 2 provides a simplification because it is the 
borderline dimension above which the rate equation ap- 
proach, dp/dt oc —p^, which gives p{t) cx 1/t, is qualita- 
tively correct because density fluctuations can be ignored 
[0, |I2|. Eq.(||) is readily adapted to this case. As be- 
fore, we treat the A particle as if it does not interact with 
the B particles, while the interactions of the B particles 
with each other give rise to their decreasing density, p{t) . 
The Poisson distribution for the number of first cross- 
ings of the A particle by B particles still holds for this 
system. The left-hand side of Eq. (|^), i.e. the probabil- 
ity density to find a B particle at the point z{t) at time 
t, becomes p(t), while on the right-hand side the propa- 
gator G{z{t),t\z{t'),t') has to be multiplied by a factor 
p{t)/p{t'), being the probability of a given B particle sur- 
viving till time t given that it survives till time t' . The 
required generalization of (nT) then reads 



' JJyA(<')G(^(^),^|^(^'),^')• 



(9) 



It is convenient to approach the limit d — > 2 from 
below. Consider first the case D' = 0, for which Q{t) 
becomes the 'site-persistence' probability, i.e. the prob- 
ability that a given point in space (the location of the 
A particle) has not been visited by any B particle. The 
static A particle corresponds to z(t) = 0, for all t, and 
with p{t') = aaiiq - l)/q]iDt')-'^^^ for large t, (|) be- 
comes 

{47rf^\a{q-l)/q= f dt' p{t')t""\t ^ t')-"'^ (10) 



for large t. In order that the right-hand side be time- 
independent for large p{t) must have the asymptotic 
form p{t) ^ Olnt. Inserting this form into (p^, and eval- 
uating the integral gives 9 = 2'^7r''/^^^ sin(7re/2)ad((7 — 
l)/q. Taking the limit e ^ 0, using = l/27re in 
this limit, gives 9 = {q — l)/q for d = 2. Finally, 
Q{t) = exp[-p{t)] ^ for t oo. The resuh 
6 = {q — l)/q in d = 2 agrees with that obtained by 
Cardy ||l^ using field-theoretic methods. 

If one now considers the case D' > 0, i.e. a diffusing 
A particle, one sees immediately that for D' = D and 
q — 2 (so that B + B ^ always), the A particle is 
equivalent to another B particle, and Q{t) = p{t), the 
density, since every surviving particle has not met any 



other particle. Hence Q{t) = p{t) \nt/t for D' ^ D 
and q = 2. This suggests that, for general D' , Q{t) will 
decay as t^^ where 6* is a nontrivial function of D' /D. 

The calculation of 9 can readily be extended to Z?' > 
within the present formalism, using a power-series expan- 
sion in D' I D. The method is to use a cumulant expan- 
sion to write Q{t) = (exp(-^))2 = cxp(-(//)^ -I- {(/x^)2 - 
(m)z}/2! + •••)■ To first order in D' , one needs only the 
first cumulant. The result \s 9 = {1 + D'/D){q - l)/q. 
This class of problems, where Q{t) decays as a power-law 
instead of a stretched exponential, marks the borderline 
where the path- integral for Q{t) is no longer dominated 
by a single path (and small fluctuations about it), giv- 
ing a leading large-t form independent of Z?', but has 
to be evaluated exactly, with results that depend on D' 
even for t ^ co. Full details of this calculation, together 
with results to higher order in D' /D, will be presented 
elsewhere. 

In conclusion, we have introduced an analytic approach 
to a class of reaction-diffusion models that reduces them 
to one-particle systems. We hope to use this method 
in future to address, inter alia, the vexing problem of 
the very slow approach to asymptopia in the trapping 
reaction. 
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